1 Detalles de la actividad


1.1 Decripción

A lo largo de esta práctica se ha realizado un caso práctico sobre el estudio de los datasets Red/White Wine Quality con el fin de aprender a identificar los datos relevantes de un proyecto analítico, así como utilizar las técnicas apropiadas para la integración, limpieza, validación y análisis de los datos.

1.2 Competencias

A continuación, se muestran las competencias desarrolladas en esta práctica del Máster de Data Science.
* Capacidad de analizar un problema en el nivel de abstracción adecuado a cada situación y aplicar las habilidades y conocimientos adquiridos para abordarlo y resolverlo.
* Capacidad para aplicar las técnicas específicas de tratamiento de datos (integración, transformación, limpieza y validación) para su posterior análisis.

1.3 Objetivos

Como se especifica en el enunciado de esta práctica, los objetivos son los siguientes:

  • Aprender a aplicar los conocimientos adquiridos y su capacidad de resolución de problemas en entornos nuevos o poco conocidos dentro de contextos más amplios o multidisciplinares.
  • Saber identificar los datos relevantes y los tratamientos necesarios (integración, limpieza y validación) para llevar a cabo un proyecto analítico.
  • Aprender a analizar los datos adecuadamente para abordar la información contenida en los datos.
  • Identificar la mejor representación de los resultados para aportar conclusiones sobre el problema planteado en el proceso analítico.
  • Actuar con los principios éticos y legales relacionados con la manipulación de datos en Tipología y ciclo de vida de los datos Práctica 2 pág 2 función del ámbito de aplicación.
  • Desarrollar las habilidades de aprendizaje que les permitan continuar estudiando de un modo que tendrá que ser en gran medida autodirigido o autónomo.
  • Desarrollar la capacidad de búsqueda, gestión y uso de información y recursos en el ámbito de la ciencia de datos.

1.4 Nota: Propiedad intelectual

A menudo es inevitable, al producir una obra multimedia, hacer uso de recursos creados por terceras personas. Es por lo tanto comprensible hacerlo en el marco de una práctica de los estudios de Informática, Multimedia y Telecomunicación de la UOC, siempre y cuando esto se documente claramente y no suponga plagio en la práctica.

Por lo tanto, si se hicera uso de recursos ajenos, se presentará un documento detallando y especificando el nombre de cada recurso, su autor, el lugar dónde se obtuvo y su estatus legal: si la obra está protegida por el copyright o se acoge a alguna otra licencia de uso (Creative Commons, licencia GNU, GPL …). Como autores de la práctica nos aseguraremos que la licencia no impide específicamente su uso en el marco de la práctica. En caso de no encontrar la información correspondiente tendrá que asumir que la obra está protegida por copyright.


2 Resolución


2.1 Introducción

2.1.1 Descripción del dataset

Como ya hemos comentado anteriormente, para el desarrollo de la práctica se ha analizado los datasets Red/White Wine Quality que podemos encontrar en UCI Machine Learning Repository. En estos ficheros cada registro representa un ensayo de calidad fisicoquímico sobre un determinado vino en base a una serie de atributos. Ambos ficheros presentan la misma configuración y el mismo número de atributos y registros como se muestra a continuación.

Formato CSV
Separador ;
Header True
Número de registros 4898
Número de atributos 12

A continuación, se muestra una descripción de los atributos de los ficheros en estudio.

Nombre Descripción Tipo
fixed acidity Acidez fija Numérico
volatile acidity Acidez volátil Numérico
citric acid Ácido cítrico Numérico
residual sugar Azúcar residual Numérico
chlorides Cloruros Numérico
free sulfur dioxide Dióxido de azufre libre Numérico
total sulfur dioxide Dióxido de azufre total Numérico
density Densidad Numérico
pH PH Numérico
sulphates Sulfatos Numérico
alcohol Porcentaje de alcohol en volumen Numérico
quality Valor de calidad resultante del ensayo Numérico

El primer paso antes de realizar el análisis será unir ambos datasets para crear un único conjunto de datos.

2.1.2 Importancia y objetivos del análisis

El objetivo es realizar un análisis completo del dataset de modo que trabajemos todas las fases de la analítica avanzada.

Tipos de análisis

Tipos de análisis

  • Análisis descriptivo. Nos permite responder a la pregunta ¿Qué ha pasado? En otras palabras, nos ofrece un resumen de los datos históricos.
  • Análisis de diagnóstico. Nos permite responder a la pregunta ¿Por qué ha pasado? En esta fase se identifican las relaciones entre las variables en estudio.
  • Análisis predictivo. Nos permite responder a la pregunta ¿Qué pasará? En esta fase se generan los modelos y se validan.
  • Análisis prescriptivo. Nos permite responder a la pregunta ¿Cómo podemos hacer que pase? Se aplica el modelo a una situación real para buscar puntos óptimos y responder a preguntas estratégicas para el diseño de nuevos productos.

2.2 Análisis del dataset

2.2.1 Instalación de librerías necesarias para exportar

library(readr) # Providing a fast and friendly way to read rectangular data ('csv'). 
library(knitr) # A General-Purpose Package for Dynamic Report Generation in R.
library(tinytex) # Helper Functions to Install and Maintain 'TeX Live', and Compile 'LaTeX' Documents. 
library(kableExtra) # Construct Complex Table with 'kable' and Pipe Syntax. Build complex HTML or 'LaTeX' tables 
library(tidyverse) # Collection of R packages designed for data science.
library(scales) # Graphical scales map data to aesthetics, and provide methods for automatically determining breaks and labels for axes
library(dplyr) # tool for working with data frame like objects, both in memory and out of memory.
library(tables) # Computes and displays complex tables of summary statistics
library(ggplot2) # A system for 'declaratively' creating graphics, based on "The Grammar of Graphics"
library(ggExtra) # Package for adding marginal histograms to ggplot2
library(magrittr) # Provides a mechanism for chaining commands with a new forward-pipe operator, %>%. 
library(corrplot) # The corrplot package is a graphical display of a correlation matrix, confidence interval....
library(grid) # grid adds an nx by ny rectangular grid to an existing plot.
library(gridExtra) # Provides a number of user-level functions to work with "grid" graphics, notably to arrange multiple grid-based plots on a page, and draw tables. 
library(visdat) # Create preliminary exploratory data visualisations of an entire dataset to identify problems or unexpected features using 'ggplot2'.
library(moments) # Compute skewness of a univariate distribution.
library(caret) # Pre-processing transformation (centering, scaling etc.) can be estimated from the training data and applied to any data set with the same variables.
library(nortest) # Performs the Lilliefors (Kolmogorov-Smirnov) test for the composite hypothesis of normality,
library(car) # VIF Regression: A Fast Regression Algorithm For Large Data
library(e1071) #Functions for latent class analysis, short time Fourier transform, fuzzy clustering, support vector machines, shortest path computation, bagged clustering, naive Bayes classifier.
library(party)# A computational toolbox for recursive partitioning.
library(rpart) # Recursive partitioning for classification, regression and survival trees
library(rpart.plot) # Plot an rpart model
library(ROCR) #ROC graphs, sensitivity/specificity curves, lift charts, and precision/recall plots
library(randomForest) #Classification and regression based on a forest of trees using random inputs, based on Breiman (2001)
library(cluster) # Methods for Cluster analysis.

2.2.2 Análisis Descriptivo

En primer lugar vamos a proceder a cargar ambos archivos. En este caso se han extraído directaemnte de la url, de modo que no tengan que pasar por nuestro dispositivo, optimizando el proceso. No obstante, se adjunta la línea de código comentada para cargar los archivos desde la carpeta raíz, por si hubiera alguna modificación de la url.

# Se carga el primer conjunto de datos
white_url <- "http://archive.ics.uci.edu/ml/machine-learning-databases/wine-quality/winequality-white.csv"
whitewine <- read.csv(white_url, header = TRUE, sep = ";")
# whitewine <- read.csv(winequality-white.csv, header = TRUE, sep = ";")

# Se carga el segundo conjunto de datos
red_url <- "http://archive.ics.uci.edu/ml/machine-learning-databases/wine-quality/winequality-red.csv"
redwine <- read.csv(white_url, header = TRUE, sep = ";")
#redwine <- read.csv(winequality-red.csv, header = TRUE, sep = ";")

#Añadimos una nueva columna con el color del vino antes de unir los dos conjuntos de datos.
whitewine$color<-"white"
redwine$color<-"red"

Unión de los conjuntos de datos

# Mediante la función rbind se unen (join) verticalmente los conjuntos de datos
wine <- rbind(redwine, whitewine)
#Una vez unidos creamos una nueva columna para futuros modelos que evaluarán la calidad
quality.factor <- ifelse(wine$quality >= 6, "Buena", "Mala")
wine <- data.frame(wine, quality.factor)
table(wine$quality.factor)
## 
## Buena  Mala 
##  6516  3280

Empezaremos haciendo un breve análisis de los datos ya que nos interesa tener una idea general de los datos que disponemos. Por ello, primero calcularemos las dimensiones de nuestra base de datos y analizaremos qué tipos de atributos tenemos.

Para empezar, calculamos las dimensiones mediante la función dim(). Obtenemos que disponemos de 9796 registros (filas) y 13 variables (columnas).

# Dimensiones de la base de datos
dim(wine)
## [1] 9796   14
# Estructura de los datos
str(wine)
## 'data.frame':    9796 obs. of  14 variables:
##  $ fixed.acidity       : num  7 6.3 8.1 7.2 7.2 8.1 6.2 7 6.3 8.1 ...
##  $ volatile.acidity    : num  0.27 0.3 0.28 0.23 0.23 0.28 0.32 0.27 0.3 0.22 ...
##  $ citric.acid         : num  0.36 0.34 0.4 0.32 0.32 0.4 0.16 0.36 0.34 0.43 ...
##  $ residual.sugar      : num  20.7 1.6 6.9 8.5 8.5 6.9 7 20.7 1.6 1.5 ...
##  $ chlorides           : num  0.045 0.049 0.05 0.058 0.058 0.05 0.045 0.045 0.049 0.044 ...
##  $ free.sulfur.dioxide : num  45 14 30 47 47 30 30 45 14 28 ...
##  $ total.sulfur.dioxide: num  170 132 97 186 186 97 136 170 132 129 ...
##  $ density             : num  1.001 0.994 0.995 0.996 0.996 ...
##  $ pH                  : num  3 3.3 3.26 3.19 3.19 3.26 3.18 3 3.3 3.22 ...
##  $ sulphates           : num  0.45 0.49 0.44 0.4 0.4 0.44 0.47 0.45 0.49 0.45 ...
##  $ alcohol             : num  8.8 9.5 10.1 9.9 9.9 10.1 9.6 8.8 9.5 11 ...
##  $ quality             : int  6 6 6 6 6 6 6 6 6 6 ...
##  $ color               : chr  "red" "red" "red" "red" ...
##  $ quality.factor      : Factor w/ 2 levels "Buena","Mala": 1 1 1 1 1 1 1 1 1 1 ...

Una vez cargados los datos, vamos a verificar la transmisión de la información.

# Hacemos una primera visualización
head(wine)
##   fixed.acidity volatile.acidity citric.acid residual.sugar chlorides
## 1           7.0             0.27        0.36           20.7     0.045
## 2           6.3             0.30        0.34            1.6     0.049
## 3           8.1             0.28        0.40            6.9     0.050
## 4           7.2             0.23        0.32            8.5     0.058
## 5           7.2             0.23        0.32            8.5     0.058
## 6           8.1             0.28        0.40            6.9     0.050
##   free.sulfur.dioxide total.sulfur.dioxide density   pH sulphates alcohol
## 1                  45                  170  1.0010 3.00      0.45     8.8
## 2                  14                  132  0.9940 3.30      0.49     9.5
## 3                  30                   97  0.9951 3.26      0.44    10.1
## 4                  47                  186  0.9956 3.19      0.40     9.9
## 5                  47                  186  0.9956 3.19      0.40     9.9
## 6                  30                   97  0.9951 3.26      0.44    10.1
##   quality color quality.factor
## 1       6   red          Buena
## 2       6   red          Buena
## 3       6   red          Buena
## 4       6   red          Buena
## 5       6   red          Buena
## 6       6   red          Buena
# Comprobamos el tipo de variable en los datos
sapply(wine,class)
##        fixed.acidity     volatile.acidity          citric.acid 
##            "numeric"            "numeric"            "numeric" 
##       residual.sugar            chlorides  free.sulfur.dioxide 
##            "numeric"            "numeric"            "numeric" 
## total.sulfur.dioxide              density                   pH 
##            "numeric"            "numeric"            "numeric" 
##            sulphates              alcohol              quality 
##            "numeric"            "numeric"            "integer" 
##                color       quality.factor 
##          "character"             "factor"
# Mostramos un resumen de los datos originales
summary(wine)
##  fixed.acidity    volatile.acidity  citric.acid     residual.sugar  
##  Min.   : 3.800   Min.   :0.0800   Min.   :0.0000   Min.   : 0.600  
##  1st Qu.: 6.300   1st Qu.:0.2100   1st Qu.:0.2700   1st Qu.: 1.700  
##  Median : 6.800   Median :0.2600   Median :0.3200   Median : 5.200  
##  Mean   : 6.855   Mean   :0.2782   Mean   :0.3342   Mean   : 6.391  
##  3rd Qu.: 7.300   3rd Qu.:0.3200   3rd Qu.:0.3900   3rd Qu.: 9.900  
##  Max.   :14.200   Max.   :1.1000   Max.   :1.6600   Max.   :65.800  
##    chlorides       free.sulfur.dioxide total.sulfur.dioxide
##  Min.   :0.00900   Min.   :  2.00      Min.   :  9.0       
##  1st Qu.:0.03600   1st Qu.: 23.00      1st Qu.:108.0       
##  Median :0.04300   Median : 34.00      Median :134.0       
##  Mean   :0.04577   Mean   : 35.31      Mean   :138.4       
##  3rd Qu.:0.05000   3rd Qu.: 46.00      3rd Qu.:167.0       
##  Max.   :0.34600   Max.   :289.00      Max.   :440.0       
##     density             pH          sulphates         alcohol     
##  Min.   :0.9871   Min.   :2.720   Min.   :0.2200   Min.   : 8.00  
##  1st Qu.:0.9917   1st Qu.:3.090   1st Qu.:0.4100   1st Qu.: 9.50  
##  Median :0.9937   Median :3.180   Median :0.4700   Median :10.40  
##  Mean   :0.9940   Mean   :3.188   Mean   :0.4898   Mean   :10.51  
##  3rd Qu.:0.9961   3rd Qu.:3.280   3rd Qu.:0.5500   3rd Qu.:11.40  
##  Max.   :1.0390   Max.   :3.820   Max.   :1.0800   Max.   :14.20  
##     quality         color           quality.factor
##  Min.   :3.000   Length:9796        Buena:6516    
##  1st Qu.:5.000   Class :character   Mala :3280    
##  Median :6.000   Mode  :character                 
##  Mean   :5.878                                    
##  3rd Qu.:6.000                                    
##  Max.   :9.000

Factorizamos la variable “color” para futuros análisis.

wine$color<-as.factor(wine$color)
wine$quality.factor<-as.factor(wine$quality.factor)

Se eliminan los registros duplicados si los hubiese

wine <- wine[!duplicated(wine), ]
dim(wine)
## [1] 7922   14

Comprobamos sie existen valores NA

vis_miss(wine)

Todos los valores están presentes. No existen valores NA.

Como podemos ver, todas las variables han sido identificadas con el tipo correspondiente y no se aprecian anomalías a simple vista en los datos.

2.2.3 Análisis Descriptivo visual

attach(wine)
#Iniciamos el cuadrante para el gráfico
par(mfrow=c(1,2))

# Representamos histogramas y boxplot para las variables numéricas    
hist(fixed.acidity , main="fixed.acidity histogram")
boxplot(fixed.acidity, main="fixed.acidity boxplot")

hist(volatile.acidity , main="volatile.acidity histogram")
boxplot(volatile.acidity , main="volatile.acidity boxplot")

hist(citric.acid  , main="citric.acid histogram")
boxplot(citric.acid  , main="citric.acid boxplot")

hist(residual.sugar  , main="residual.sugar histogram")
boxplot(residual.sugar  , main="residual.sugar boxplot")

hist(chlorides  , main="chlorides histogram")
boxplot(chlorides  , main="chlorides boxplot")

hist(free.sulfur.dioxide  , main="free.sulfur.dioxide histogram")
boxplot(free.sulfur.dioxide  , main="free.sulfur.dioxide boxplot")

hist(total.sulfur.dioxide  , main="total.sulfur.dioxide histogram")
boxplot(total.sulfur.dioxide  , main="total.sulfur.dioxide boxplot")

hist(density  , main="density histogram")
boxplot(density  , main="density boxplot")

hist(sulphates  , main="sulphates histogram")
boxplot(sulphates  , main="sulphates boxplot")

hist(alcohol  , main="alcohol histogram")
boxplot(alcohol  , main="alcohol boxplot")

hist(quality  , main="quality histogram")
boxplot(quality  , main="quality boxplot")

# Representamos gráfico de barras de las variables nominales
barplot(table(color), main="Color histogram")
barplot(table(quality.factor), main="quality.factor histogram")

detach(wine)

Todas las variables presentan outliers.

  • Quality tiene la mayoría de los valores concentrados en las categorías 5, 6 y 7. Sólo una pequeña proporción está en las categorías [3, 4] y [8, 9] y ninguna en las categorías [1, 2, 10].
  • Fixed acidity, volatile.acidity y citric.acid tienen valores atípicos. Si se eliminan esos valores atípicos, la distribución de las variables puede considerarse simétrica.
  • Residual.sugar residual tiene una distribución positivamente sesgada; incluso después de eliminar los valores atípicos, la distribución seguirá siendo sesgada.
  • Algunas de las variables, por ejemplo, free.sulphur.dioxide y density,, tienen unos pocos valores atípicos, pero éstos son muy diferentes del resto.
  • Alcohol tiene una distribución irregular pero no tiene valores atípicos pronunciados.

2.2.4 Comprobación de la asimetría

Comprobamos a continuación la asimetría de las variables verificar si los datos se distribuyen normalmente o no.

attach(wine)
skewness(fixed.acidity)
## [1] 0.6957048
skewness(volatile.acidity)
## [1] 1.640149
skewness(citric.acid)
## [1] 1.309857
skewness(residual.sugar)
## [1] 1.332882
skewness(chlorides)
## [1] 4.966254
skewness(free.sulfur.dioxide)
## [1] 1.56579
skewness(total.sulfur.dioxide)
## [1] 0.4565402
skewness(density)
## [1] 1.272595
skewness(pH)
## [1] 0.4551981
skewness(sulphates)
## [1] 0.9373206
skewness(alcohol)
## [1] 0.4504406
skewness(quality)
## [1] 0.1119404
detach(wine)

Regla de la asimetría:

  • Si la asimetría es 0, los datos son perfectamente simétricos.
  • Si la asimetría es menor que -1 o mayor que 1, la distribución es altamente sesgada.
  • Si la asimetría está entre -0,5 y 0,5, la distribución es aproximadamente simétrica.

Fuente:*GoodData.Documentation

2.2.5 Transformación de las variables

Se aplica transformación de Boxcox y volvemos a comprobar la asimetría.

Fuente: *Box-Cox Transformations

wine_procesado <- preProcess(wine[,1:11], c("BoxCox", "center", "scale"))
new_wine <- data.frame(trans = predict(wine_procesado, wine))

colnames(new_wine)
##  [1] "trans.fixed.acidity"        "trans.volatile.acidity"    
##  [3] "trans.citric.acid"          "trans.residual.sugar"      
##  [5] "trans.chlorides"            "trans.free.sulfur.dioxide" 
##  [7] "trans.total.sulfur.dioxide" "trans.density"             
##  [9] "trans.pH"                   "trans.sulphates"           
## [11] "trans.alcohol"              "trans.quality"             
## [13] "trans.color"                "trans.quality.factor"
attach(new_wine)
skewness(trans.fixed.acidity)
## [1] 0.1131602
skewness(trans.volatile.acidity)
## [1] 0.1718966
skewness(trans.citric.acid)
## [1] 1.309857
skewness(trans.residual.sugar)
## [1] -0.05827602
skewness(trans.chlorides)
## [1] -0.1496442
skewness(trans.free.sulfur.dioxide)
## [1] 0.0950526
skewness(trans.total.sulfur.dioxide)
## [1] 0.008932695
skewness(trans.density)
## [1] 1.099868
skewness(trans.quality)
## [1] 0.1119404
skewness(trans.pH)
## [1] 0.0002881811
skewness(trans.sulphates)
## [1] -0.01873644
skewness(trans.alcohol)
## [1] 0.03929748
detach(new_wine)

Si el coeficiente de asimetría no se aproxima a 0, entonces su conjunto de datos no se distribuye normalmente.

Ahora la mayoría las distribuciones son aproximadamente simétricas (la asimetría está entre -0,5 y 0,5). Las únicas variables que muestran aún una asimestría son:
* trans.citric.acid * trans.density

Se a va a comprobr mediante gráfico de cuantiles teóricos (Gráficos Q-Q) la normalidad de las dos variables que muestran aín coeficientes de simetría alejados de 0.
Para el resto de variables no es necesario comprobarlo ya que sus coeficientes se distribuyen alrededor de 0 y por tanto se distribuyen normalmente.

Se realiza una comprobación de la normalidad para la variable trans.citric.acid.

# Para comprobar la normalidad de la variable trans.citric.acid se va a usar:
# 1. Gráfico de cuantiles teóricos (Gráficos Q-Q): Consiste en comparar los cuantiles de la distribución observada con los cuantiles teóricos de una distribución normal con la misma media y desviación estándar que los datos. Cuanto más se aproximen los datos a una normal, más alineados están los puntos entorno a la recta.
# 2. Un histograma
# 3. Diagrama de cajas
par(mfrow = c(1, 3)) 
hist(new_wine$trans.citric.acid, las=1, main="Distribución normal", font.main=1, ylab="Frecuencia", xlab = "trans.citric.acid") 
qqnorm(new_wine$trans.citric.acid, las=1, pch=18, main="Simetria", font.main=1, xlab="Cuantiles teóricos", ylab="Cuantiles muestrales") 
qqline(new_wine$trans.citric.acid) 
boxplot(new_wine$trans.citric.acid, las=1, main="Valores extremos", font.main=1, xlab="Datos de distribution normal", horizontal=F)

No se trata de una distribución normal ya que presenta una asimetría en ambas colas, aunque la cola a la derecha de la media es más larga que la izquierda.

Para asegurarnos completamente se realiza la pruena de Lilliefors (Kolmogorov-Smirnov) para muestras grandes (n>50).
El test Lilliefors asume que la media y varianza son desconocidas, estando especialmente desarrollado para contrastar la normalidad.
Es la alternativa al test de Shapiro-Wilk cuando el número de observaciones es mayor de 50. La función lillie.test() del paquete nortest permite aplicarlo.

Planteamos la Hipótesis:

\(\ H_0:\) La muestra proviene de una distribución normal.
\(\ H_1:\) La muestra no proviene de una distribución normal.

Nivel de significancia: El nivel de significancia que se trabajará es de 0.05. \(\alpha=0.05\)

Criterio de decisón:

Si \(p<\alpha\) se rechaza \(\ H_0:\)
Si \(p>\alpha\) no se rechaza \(\ H_0:\)

lillie.test(x=new_wine$trans.citric.acid)
## 
##  Lilliefors (Kolmogorov-Smirnov) normality test
## 
## data:  new_wine$trans.citric.acid
## D = 0.10992, p-value < 2.2e-16

Dado que \(p-value<\alpha\) se rechaza \(\ H_0:\). La muestra no proviene de una distribución normal. Pero, puesto que el tamaño de la muestra es mayor a 30, según el teorema del límite central, asumimos la normalidad.

Se realiza una comprobación de la normalidad para la variable trans.density.

# Para comprobar la normalidad de la variable trans.citric.acid se va a usar:
# 1. Gráfico de cuantiles teóricos (Gráficos Q-Q): Consiste en comparar los cuantiles de la distribución observada con los cuantiles teóricos de una distribución normal con la misma media y desviación estándar que los datos. Cuanto más se aproximen los datos a una normal, más alineados están los puntos entorno a la recta.
# 2. Un histograma
# 3. Diagrama de cajas
par(mfrow = c(1, 3)) 
hist(new_wine$trans.density, las=1, main="Distribución normal", font.main=1, ylab="Frecuencia", xlab = "trans.density") 
qqnorm(new_wine$trans.density, las=1, pch=18, main="Simetria", font.main=1, xlab="Cuantiles teóricos", ylab="Cuantiles muestrales") 
qqline(new_wine$trans.density) 
boxplot(new_wine$trans.density, las=1, main="Valores extremos", font.main=1, xlab="Datos de distribution normal", horizontal=F)

Presenta una asimetría en en la cola izquierda.

Para asegurarnos completamente se realiza la pruena de Lilliefors (Kolmogorov-Smirnov) para muestras grandes (n>50).
El test Lilliefors asume que la media y varianza son desconocidas, estando especialmente desarrollado para contrastar la normalidad.
Es la alternativa al test de Shapiro-Wilk cuando el número de observaciones es mayor de 50. La función lillie.test() del paquete nortest permite aplicarlo.

Planteamos la Hipótesis:

\(\ H_0:\) La muestra proviene de una distribución normal.
\(\ H_1:\) La muestra no proviene de una distribución normal.

Nivel de significancia: El nivel de significancia que se trabajará es de 0.05. \(\alpha=0.05\)

Criterio de decisón:

Si \(p<\alpha\) se rechaza \(\ H_0:\)
Si \(p>\alpha\) no se rechaza \(\ H_0:\)

lillie.test(x=new_wine$trans.density)
## 
##  Lilliefors (Kolmogorov-Smirnov) normality test
## 
## data:  new_wine$trans.density
## D = 0.048519, p-value < 2.2e-16

Dado que \(p-value<\alpha\) se rechaza \(\ H_0:\). La muestra no proviene de una distribución normal. Pero, puesto que el tamaño de la muestra es mayor a 30, según el teorema del límite central, asumimos la normalidad.

2.2.6 Eliminación de los Outliers

Se analiza cada variable independientemente para determinar qué valores atípicos, si los hay, deben ser eliminados.

trans.fixed.acidity

boxplot(trans.fixed.acidity ~ trans.quality , data=new_wine,
        xlab="Quality", ylab="Fixed Acidity")

summary(new_wine$trans.fixed.acidity)
##     Min.  1st Qu.   Median     Mean  3rd Qu.     Max. 
## -4.64627 -0.59543  0.01653  0.00000  0.58505  5.91650
quantile(new_wine$trans.fixed.acidity, c(.999))
##    99.9% 
## 3.343608
new_wine$trans.quality[new_wine$trans.fixed.acidity > 3.4]
## [1] 6 6 3 6 6 3

Se observan bastantes valores atípicos, especialmente para la calidad de 6. El mayor es tres veces más grande que el valor del tercer cuartil.
El 99,9% de las observaciones tienen una acidez fija tras la normalización de los datos de 3,34 o menos.
Hay 6 outliers. De nuevo, la mayoría de ellas están clasificadas con una calidad de 6. Esto no es sorprendente, porque muchas observaciones en el conjunto de datos están clasificadas como calidad 6.

trans.volatile.acidity

boxplot(trans.volatile.acidity ~ trans.quality , data=new_wine,
        xlab="Quality", ylab="volatile acidity")

summary(new_wine$trans.volatile.acidity)
##     Min.  1st Qu.   Median     Mean  3rd Qu.     Max. 
## -3.49849 -0.67297 -0.04768  0.00000  0.65033  4.17526
quantile(new_wine$trans.volatile.acidity, c(.999))
##    99.9% 
## 3.625124
new_wine$trans.quality[new_wine$trans.volatile.acidity > 3.7]
## [1] 4 6 4 4 6 4

Esta vez la mayoría de los outliers tienen la calidad 4. Esto es interesante porque no hay muchas observaciones con la calidad 4. La acidez puede ser un factor clave para predecir los vinos malos, por ello puede que estos no sean en realidad valores atípicos.

trans.citric.acid

boxplot(trans.citric.acid ~ trans.quality , data=new_wine,
        xlab="Quality", ylab="citric acid")

summary(new_wine$trans.citric.acid)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
## -2.7306 -0.5254 -0.1171  0.0000  0.4547 10.8272
quantile(new_wine$trans.citric.acid, c(.999))
##    99.9% 
## 5.436758
new_wine$trans.quality[new_wine$trans.citric.acid > 5.5]
## [1] 6 6 6 6

La mayoría de outliers están clasificados como 6 en calidad. El boxplot es interesante porque tanto el vino de calidad 3 como el de calidad 9 no tienen outliers y todos los demás rangos sí. Esto podría ser sólo una casualidad. También podría ser que la variación en el ácido cítrico está de alguna manera relacionada con los vinos de nivel medio.

trans.residual.sugar

boxplot(trans.residual.sugar ~ trans.quality , data=new_wine,
        xlab="Quality", ylab="residual sugar")

summary(new_wine$trans.residual.sugar)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
## -2.1014 -1.0237  0.1603  0.0000  0.8619  3.0601
quantile(new_wine$trans.residual.sugar, c(.999))
##    99.9% 
## 1.889274
new_wine$trans.quality[new_wine$trans.residual.sugar > 2.0]
## [1] 6 6 6 6 6 6

EL boxplot no muestra valores atípicos. La calidad 6 muestra 6 valores altos.

trans.chlorides

boxplot(trans.chlorides ~ trans.quality , data=new_wine,
        xlab="Quality", ylab="chlorides")

summary(new_wine$trans.chlorides)
##     Min.  1st Qu.   Median     Mean  3rd Qu.     Max. 
## -7.19158 -0.56676  0.02707  0.00000  0.54648  4.08105
quantile(new_wine$trans.chlorides, c(.999))
##    99.9% 
## 3.729705
new_wine$trans.quality[new_wine$trans.chlorides > 3.8]
## [1] 5 4 5 5 4 5

Los vinos con calidad 9 están muy bien distribuidos y no hay outliers. Además, la media es más baja para este grupo que para cualquier otro. Es posible que un nivel muy bajo de cloruros sea indicativo de un muy buen vino. También podría ser debido al número extremadamente pequeño de vinos calidad 9.

trans.free.sulfur.dioxide

boxplot(trans.free.sulfur.dioxide ~ trans.quality , data=new_wine,
        xlab="Quality", ylab="free sulfur dioxide")

summary(new_wine$trans.free.sulfur.dioxide)
##     Min.  1st Qu.   Median     Mean  3rd Qu.     Max. 
## -2.96992 -0.64046  0.01309  0.00000  0.67690  7.76653
quantile(new_wine$trans.free.sulfur.dioxide, c(.999))
##    99.9% 
## 3.856631
new_wine$trans.quality[new_wine$trans.free.sulfur.dioxide > 3.9]
## [1] 5 3 4 3 5 3 4 3

Todos las calidades tienen al menos un atípico y los medios son bastante similares con la excepción de 4. Parece que a medida que la calidad mejora, el rango de variación tiende a disminuir.

trans.total.sulfur.dioxide

boxplot(trans.total.sulfur.dioxide ~ trans.quality , data=new_wine,
        xlab="Quality", ylab="total sulfur dioxide")

summary(new_wine$trans.total.sulfur.dioxide)
##     Min.  1st Qu.   Median     Mean  3rd Qu.     Max. 
## -4.17206 -0.69046 -0.03296  0.00000  0.70174  5.39107
quantile(new_wine$trans.total.sulfur.dioxide, c(.999))
##    99.9% 
## 3.349804
new_wine$trans.quality[new_wine$trans.total.sulfur.dioxide > 3.5]
## [1] 3 5 3 3 5 3

La variación tiende a disminuir con la mejora de la calidad del vino.

trans.density

boxplot(trans.density ~ trans.quality , data=new_wine,
        xlab="Quality", ylab="trans density")

summary(new_wine$trans.density)
##     Min.  1st Qu.   Median     Mean  3rd Qu.     Max. 
## -2.33088 -0.74894 -0.09586  0.00000  0.66713 14.63895
quantile(new_wine$trans.density, c(.999))
##    99.9% 
## 2.946401
new_wine$trans.quality[new_wine$trans.density > 3.0]
## [1] 6 6 6 6 6 6

La calidad 6 muestra un claro outlier. Dada la cantidad de vinos con calidad 6, parece poco probable que una densidad de este nivel sea precisa o representativa del grupo.

trans.pH

boxplot(trans.pH ~ trans.quality , data=new_wine,
        xlab="Quality", ylab="trans pH")

summary(new_wine$trans.pH)
##     Min.  1st Qu.   Median     Mean  3rd Qu.     Max. 
## -3.95467 -0.67121 -0.03552  0.00000  0.67393  3.30014
quantile(new_wine$trans.pH, c(.999))
##    99.9% 
## 3.182558
new_wine$trans.quality[new_wine$trans.pH > 3.2]
## [1] 7 6 6 6 7 6 6 6

Se observa una mayor dispersión en la calidad 3, pero no hay valores atípicos, mientras que la calidad 9 muestra una menor dispersión.

trans.sulphates

boxplot(trans.sulphates ~ trans.quality , data=new_wine,
        xlab="Quality", ylab="trans sulphates")

summary(new_wine$trans.sulphates)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
## -3.8891 -0.6722  0.0514  0.0000  0.6494  3.2783
quantile(new_wine$trans.sulphates, c(.999))
##    99.9% 
## 2.971275
new_wine$trans.quality[new_wine$trans.sulphates > 3.0]
## [1] 6 6 6 7 6 6 6 7

Las medias de las calidades estaán todas en el mismo nivel. Es difícil diferenciar, a partir de las medias, si los sulfatos pueden llegar a influenciar en la calidad.

trans.alcohol

boxplot(trans.alcohol ~ trans.quality , data=new_wine,
        xlab="Quality", ylab="trans alcohol")

summary(new_wine$trans.alcohol)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
## -2.8604 -0.8998 -0.0304  0.0000  0.7529  2.2877
quantile(new_wine$trans.alcohol, c(.999))
##    99.9% 
## 2.201294
new_wine$trans.quality[new_wine$trans.alcohol > 2.21]
## [1] 7 7 7 7
new_wine$trans.quality[new_wine$trans.alcohol < -2.8]
## [1] 5 3 5 3

Parece haber una relación entre el alcohol y la calidad.

Procedemos a eliminar los valores outliers encontardos en los gráficos boxplot anteriores.

# Mostramos los valores identificados como outliers
attach(new_wine)
f_a_out<-boxplot.stats(trans.fixed.acidity)$out
v_a_out<-boxplot.stats(trans.volatile.acidity)$out
c_a_out<-boxplot.stats(trans.citric.acid)$out
r_s_out<-boxplot.stats(trans.residual.sugar)$out
ch_out<-boxplot.stats(trans.chlorides)$out
f_s_d_out<-boxplot.stats(trans.free.sulfur.dioxide)$out
t_s_d_out<-boxplot.stats(trans.total.sulfur.dioxide)$out
d_t_out<-boxplot.stats(trans.density)$out
ph_out<-boxplot.stats(trans.pH )$out
sh_out<-boxplot.stats(trans.sulphates)$out
ac_out<-boxplot.stats(trans.alcohol)$out
qa_out<-boxplot.stats(trans.quality)$out
detach(new_wine)

Las variables trans.alcohol, trans.residual.sugar y trans.color no presentan outliers.

# Eliminamos las instancias con outlier
new_wine<-new_wine[-which(new_wine$trans.fixed.acidity %in% f_a_out),]
new_wine<-new_wine[-which(new_wine$trans.volatile.acidity %in% v_a_out),]
new_wine<-new_wine[-which(new_wine$trans.citric.acid %in% c_a_out),]
new_wine<-new_wine[-which(new_wine$trans.chlorides %in% ch_out),]
new_wine<-new_wine[-which(new_wine$trans.free.sulfur.dioxide %in% f_s_d_out),]
new_wine<-new_wine[-which(new_wine$trans.total.sulfur.dioxide %in% t_s_d_out),]
new_wine<-new_wine[-which(new_wine$trans.density %in% d_t_out),]
new_wine<-new_wine[-which(new_wine$trans.pH %in% ph_out),]
new_wine<-new_wine[-which(new_wine$trans.sulphates %in% sh_out),]
new_wine<-new_wine[-which(new_wine$trans.quality %in% qa_out),]

2.2.7 Análisis de diagnóstico

Tras haber eliminado los valores atípicos, revisamos la correlación entre las distintas variables. El objetivo es centrar el análisis descriptivo en aquellas variables que presentan correlación.

#Se omite la variable color para pocer representar la matriz de correlación
corrplot(cor(new_wine[1:12]))

Las correlaciones positivas se muestran en azul y las negativas en rojo. La intensidad del color y el tamaño del círculo son proporcionales a los coeficientes de correlación. En la parte derecha del correlograma, el color de la leyenda muestra los coeficientes de correlación y los colores correspondientes.

Observando los resultados se verifica que:
1. La calidad del vino no se encuentra correlacionada con las variables “citric.acid”, “free.sulfur.dioxide”, los “sulphates” y “color”. 2. Es coherente la alta correlación positiva entre free.sulfur.dioxide y total.sulfur.dioxide. A mayor cantidad de sulfuro más sulfuro libre hay. 3. Es coherente la alta correlación negativa entre fixed.acidity y PH. Cuanto mayor es el PH en la escala de ácido-base más básico es el compuesto.Por lo tanto a valores altos de PH menor es la acidez. 4. Es coherente la alta correlación negativa entre el azucar y el alcohol.El alcohol es el resultado de la oxidación del azucar, por tanto si hay menos azucar el grado de alcohol aumenta.

A continuación se estudian visaulamente las posibles correlaciones entre:

2.2.8 Variables residual.sugar y density

# Calculamos la correlación y la mostramos por pantalla
RS_D_cor = cor(new_wine$trans.residual.sugar,new_wine$trans.density)
print(RS_D_cor)
## [1] 0.76787
# Generamos el scatter plot y añadimos la linea de regresión para una mejor visualización
plot(new_wine$trans.residual.sugar,new_wine$trans.density, main="residual.sugar vs density",xlab="residual.sugar", ylab="density", pch=19)
abline(lm(new_wine$trans.density~new_wine$trans.residual.sugar), col="red") # regression line (y~x)

Como podemos apreciar en la gráfica existe una fuerte relación lineal entre ambas variables y todos los puntos se ajustan bastante bien a la linea. El resultado de la correlación, 0.75, refuerza lo comentado anteriormente.

2.2.9 Variables density y alcochol

# Calculamos la correlación y la mostramos por pantalla
D_A_cor = cor(new_wine$trans.density,new_wine$trans.alcohol)
print(D_A_cor)
## [1] -0.806588
# Generamos el scatter plot y añadimos la linea de regresión para una mejor visualización
plot(new_wine$trans.density,new_wine$trans.alcohol, main="alcochol vs density",xlab="alcochol", ylab="density", pch=19)
abline(lm(new_wine$trans.alcohol~new_wine$trans.density), col="red") # regression line (y~x)

Al igual que antes, existe una fuerte relación lineal entre ambas variables y los puntos se ajustan a la linea. El resultado de la correlación, -0.80, refuerza lo comentado anteriormente.

2.2.10 Todas las variables y quality

# Calculamos la correlación y la mostramos por pantalla
par(mfrow = c(1,2))
for (i in c(1:11)) {
    plot(new_wine[, i], jitter(new_wine[, "trans.quality"]), xlab = names(new_wine)[i],
         ylab = "trans.quality", col = "firebrick", cex = 0.8, cex.lab = 1.3)
    abline(lm(new_wine[, "trans.quality"] ~ new_wine[ ,i]), lty = 2, lwd = 2)
}

par(mfrow = c(1, 1))

Las variables trans.alcohol, trans.sulphates y trans.pH parecen tener una correlación positiva con trans.quality, mientras que ltrans.volatile.acidity, trans.chlorides y trans.total.sulfur.dioxidel parecen tener una relación negativa con la calidad.

2.2.11 Contraste de hipótesis para la diferencia de medias

Dado que nuestro dataset contiene los valores para los vinos blancos y tintos, vamos a realizar un contraste de hipótesis para la diferencia de medias (trans.citric.acid) de las dos muestras. Planteamos la sigueinte hipótesis:
¿Podemos aceptar con un nivel de confianza del 95% que los vinos blancos (color) tienen una acidez (trans.citric.acid) que supera a la acidez de los vinos tintos (color)?.

Planteamos la Hipótesis nula y alternativa \[ \left\{ \begin{array}{ll} H_{0}: & \mu_{si}=\mu_{no}\\ H_{1}: & \mu_{si}>\mu_{no} \end{array} \right.\\ \ Hipótesis\ unilateral\ a \ la \ dereccha \]

EL método que se aplica es un contraste de hipótesis para dos muestras independientes. El test es paramétrico, asumiendo el teorema del límite central. Es un test unilateral, con varianzas desconocidas, pero iguales, ya que el test de varianzas no rechaza la hipótesis nula de varianzas iguales.

#Test de igualdad de varianzas
#H0: F (ratio de varianzas) = 1
#H1: F diferente de 1
var.test(new_wine$trans.citric.acid[new_wine$trans.color=="white"], new_wine$trans.citric.acid[new_wine$trans.color=="red"],
         alternative = "greater")
## 
##  F test to compare two variances
## 
## data:  new_wine$trans.citric.acid[new_wine$trans.color == "white"] and new_wine$trans.citric.acid[new_wine$trans.color == "red"]
## F = 1, num df = 3242, denom df = 3242, p-value = 0.5
## alternative hypothesis: true ratio of variances is greater than 1
## 95 percent confidence interval:
##  0.9438529       Inf
## sample estimates:
## ratio of variances 
##                  1

A continuación se lleva a cabo el contraste de hipótesis:

#Se define una función myttest.
#x1, x2: muestras
#CL: confidence level (nivel de confianza)
#equalvar TRUE/FALSE
#alternative: "bilateral", "less" (x1 less than x2), "greater" (x1 greater than x2).
myttest <- function( x1, x2, CL=0.95,equalvar=TRUE, alternative="bilateral" ){ # z test
  mean1<-mean(x1)
  n1<-length(x1)
  sd1<-sd(x1)
  mean2<-mean(x2) 
  n2<-length(x2) 
  sd2<-sd(x2)
  
if (equalvar==TRUE){
  s <-sqrt( ((n1-1)*sd1^2 + (n2-1)*sd2^2 )/(n1+n2-2) ) 
  Sb <- s*sqrt(1/n1 + 1/n2)
  df<-n1+n2-2
}
else{ #equalvar==FALSE
  Sb <- sqrt( sd1^2/n1 + sd2^2/n2 )
  denom <- ( (sd1^2/n1)^2/(n1-1) + (sd2^2/n2)^2/(n2-2)) 
  df <- ((sd1^2/n1 + sd2^2/n2)^2) / denom
}
alfa <- (1-CL)
t<- (mean1-mean2) / Sb

if (alternative=="bilateral"){
  tcritical <- qt( alfa/2, df, lower.tail=FALSE ) #two sided 
  pvalue<-pt( abs(t), df, lower.tail=FALSE )*2 #two sided
}
else if (alternative=="less"){
  tcritical <- qt( alfa, df, lower.tail=TRUE )
  pvalue<-pt( t, df, lower.tail=TRUE ) 
}

else{ #(alternative=="greater")
  tcritical <- qt( alfa, df, lower.tail=FALSE ) 
  pvalue<-pt( t, df, lower.tail=FALSE )
}
  #Guardamos en resultado en un data frame
info<-data.frame(t,df,tcritical,pvalue) 
info %>% kable() %>% kable_styling() 
return (info)
}

Se visualizan tanto el estadístico como el valor p:

info<-myttest( new_wine$trans.citric.acid[new_wine$trans.color=="white"], 
               new_wine$trans.citric.acid[new_wine$trans.color=="red"], 
               alternative="less",
               equalvar="TRUE")

info 
##   t   df tcritical pvalue
## 1 0 6484 -1.645089    0.5

Dado que hemos aplicado una función propia comprobamos el resultado con la fución original de R.

#Comparación con t.test de R
t.test( new_wine$trans.citric.acid[new_wine$trans.color=="white"], new_wine$trans.citric.acid[new_wine$trans.color=="red"], alternative = "less",
        var.equal = TRUE)
## 
##  Two Sample t-test
## 
## data:  new_wine$trans.citric.acid[new_wine$trans.color == "white"] and new_wine$trans.citric.acid[new_wine$trans.color == "red"]
## t = 0, df = 6484, p-value = 0.5
## alternative hypothesis: true difference in means is less than 0
## 95 percent confidence interval:
##        -Inf 0.02953236
## sample estimates:
##   mean of x   mean of y 
## -0.07990945 -0.07990945

El valor p es 0.5, mayor que el nivel de significación es 0.05 (95%). No podemos rechazar la hipótesis nula de igualdad de acidez media entre los vinos blancos y tintos. En otras palabras, no podemos asumir que la media de acidez de vinos blancos sea mayor a la de vinos tintos.

2.2.12 Contraste de hipótesis para la relación entre variables

  • Hipótesis nula: H0: La variable trans.color y la variable trans.quality son independientes.
  • Hipótesis alternativa: H1: La variable trans.color y la variable trans.quality no son independientes.

Realizad los cálculos manualmente del test. Para ello, se recomienda construir una función my.chisq que realice estos cálculos y que podáis reusar más adelante.

# Realizamos el test con la función chisq.test para comprobar los resultados.
chisq.test(new_wine$trans.color,new_wine$trans.quality)
## 
##  Pearson's Chi-squared test
## 
## data:  new_wine$trans.color and new_wine$trans.quality
## X-squared = 0, df = 3, p-value = 1

El valor crítico es mayor que el chi-cuadrado, por tanto, no podemos rechazar la hipótesis nula de que las variables son independientes.

Del mismo modo, el p valor es mayor al nivel de confianza 0.05 por tanto también debemos aceptar la hipótesis nula y podemos afirmar que no hay una dependencia entre el color y si la calidad del vino con un 95% de confianza.

2.2.13 Análisis predictivo

Vamos a aplicar un modelo de regresión lineal multivariable. El objetivo consistirá en predecir la calidad del vino en función de las otras variables.

Antes de elegir el mejor modelo lineal vamos a dividir el conjunto de datos en dos:
1. Set de entrenamiento: Un subconjunto para entrenar un modelo (80% de los datos)
2. Set de prueba:Un subconjunto para probar el modelo entrenado (20% de lod datos)

División del conjunto de datos

División del conjunto de datos

Esta división se ha llevado a cabo para cumplir con los siguientes condiciones: * Que sea lo suficientemente grande como para generar resultados significativos desde el punto de vista estadístico.
* Que sea representativo de todo el conjunto de datos.

Se divide el conjunto de datos (80% entrenamiento - 20% test)

set.seed(100)
#new_wine$trans.color<-as.numeric(new_wine$trans.color)
trainingRowIndex <- sample(1:nrow(new_wine), 0.8*nrow(new_wine))  
new_wine_train <- new_wine[trainingRowIndex, ]  
new_wine_test  <- new_wine[-trainingRowIndex, ]  

Regresión lineal Múltiple.

Realizamos el primer modelo

# Incluimos todas las variables 
modelo_1 <- lm(trans.quality ~ . , new_wine_train) 
summary(modelo_1)
## 
## Call:
## lm(formula = trans.quality ~ ., data = new_wine_train)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.23450 -0.26344 -0.07044  0.16755  1.07462 
## 
## Coefficients:
##                             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                 6.232979   0.008759 711.633  < 2e-16 ***
## trans.fixed.acidity         0.043538   0.008573   5.079 3.93e-07 ***
## trans.volatile.acidity     -0.011169   0.006488  -1.721 0.085232 .  
## trans.citric.acid           0.017253   0.007862   2.195 0.028240 *  
## trans.residual.sugar        0.118971   0.013485   8.822  < 2e-16 ***
## trans.chlorides            -0.032990   0.008514  -3.875 0.000108 ***
## trans.free.sulfur.dioxide   0.059832   0.007766   7.704 1.57e-14 ***
## trans.total.sulfur.dioxide -0.020586   0.008534  -2.412 0.015894 *  
## trans.density              -0.193877   0.022612  -8.574  < 2e-16 ***
## trans.pH                    0.053138   0.007315   7.264 4.32e-13 ***
## trans.sulphates             0.027806   0.006044   4.600 4.32e-06 ***
## trans.alcohol              -0.003995   0.014352  -0.278 0.780729    
## trans.colorwhite           -0.003302   0.010646  -0.310 0.756462    
## trans.quality.factorMala   -1.241168   0.012944 -95.888  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.3834 on 5174 degrees of freedom
## Multiple R-squared:  0.7437, Adjusted R-squared:  0.7431 
## F-statistic:  1155 on 13 and 5174 DF,  p-value: < 2.2e-16

Todas las variables presentan significación muy alta, esto puede ser consecuencia de multicolinealidad entre las variables. Para detectarla y minimizarla vamos a seguir los siguientes pasos:

  1. Análisis de la multicolinealidad
  2. Modelo con variables que no muestran colinealidad
  3. Análisis de la multicolinealidad
  4. Modelo con variables más significativas
  5. Predicción de valores
  6. Matriz de confusión

Multicolinealidad: Se dice que existe multicolinealidad entre las variables explicativas cuando existe algún tipo de dependencia lineal entre ellas, o lo que es lo mismo, si existe una fuerte correlación entre las mismas. La correlación no solamente se refiere a las distintas variables dos a dos, sino a cualquier de ellas con cualquier grupo de las restantes. Por esta razón no es suficiente (aunque sí necesaria) que en la matriz de correlaciones bivariadas haya correlaciones altas.

El principal inconveniente de la multicolinealidad consiste en que se incrementan la varianza de los coeficientes de regresión estimados hasta el punto que resulta prácticamente imposible establecer su significación estadística, ya que como se sabe, el valor de t para un determinado coeficiente de regresión es el valor de dicho coeficiente dividido por su desviación tipica. Si este es grande, el valor de t será bajo y no llegara a la significación.

Hay varios procedimientos para detectar la multicolinealidad entre los predictores. El primero de ellos, basado en la correlación múltiple de un determinado regresor con los restantes se denomina Tolerancia de dicho regresor. Su valor es \(\ 1-R^2_i\). Siendo \(\ R^2_i\) la correlación múltiple al cuadrado de dicho regresor con los restantes.

Para que haya multicolinealidad dicha correlación ha de ser alta, o lo que es lo mismo la tolerancia baja.

Otro índice relacionado con éste y que nos da una idea del grado de aumento de la varianza se denomina Factor de Inflación de la Varianza, y es precisamente el recíproco de la tolerancia. Su valor es: \(\ VIF=1/(1-R^2_y)\). Para que no haya multicolinealidad el denominador tiene que valer cerca de la unidad, por tanto un poco más de 1 el valor de VIF. Cuanto mayor sea de este valor mayor multicolinealidad habrá y por tanto más inestable es el modelo lineal; en otras palabras, es conveniente eliminar del modelo las variables con un factor VIF grande.

Si todos los VIF son 1, no hay multicolinealidad, pero si algunos VIF son mayores que 1, los predictores están correlacionados. Cuando un VIF es > 5, el coeficiente de regresión para ese término no se estima adecuadamente.

Fuentes: REGRESIÓN LINEAL MÚLTIPLE
Multicolinealidad

Analizamos la multicolinealidad mediante el análisis del valor VIF

vif(modelo_1)
##        trans.fixed.acidity     trans.volatile.acidity 
##                   2.058439                   1.234414 
##          trans.citric.acid       trans.residual.sugar 
##                   1.133725                   6.336802 
##            trans.chlorides  trans.free.sulfur.dioxide 
##                   1.466964                   1.871523 
## trans.total.sulfur.dioxide              trans.density 
##                   2.375456                  16.440751 
##                   trans.pH            trans.sulphates 
##                   1.658096                   1.115524 
##              trans.alcohol                trans.color 
##                   6.767601                   1.000318 
##       trans.quality.factor 
##                   1.305766

Se observa multicolinealidad en la variable trans.density y se procede a eliminarla del modelo.
Ajustamos el modelo de nuevo

# excluímos la variable trans.density
modelo_2 <- lm(trans.quality ~ . -trans.density, new_wine_train) 
summary(modelo_2)
## 
## Call:
## lm(formula = trans.quality ~ . - trans.density, data = new_wine_train)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.20587 -0.26420 -0.08301  0.16842  1.07027 
## 
## Coefficients:
##                             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                 6.236819   0.008808 708.061  < 2e-16 ***
## trans.fixed.acidity        -0.001727   0.006801  -0.254 0.799617    
## trans.volatile.acidity     -0.007167   0.006517  -1.100 0.271503    
## trans.citric.acid           0.011003   0.007882   1.396 0.162790    
## trans.residual.sugar        0.017403   0.006489   2.682 0.007345 ** 
## trans.chlorides            -0.042021   0.008508  -4.939 8.09e-07 ***
## trans.free.sulfur.dioxide   0.063904   0.007806   8.187 3.34e-16 ***
## trans.total.sulfur.dioxide -0.027176   0.008559  -3.175 0.001506 ** 
## trans.pH                    0.021751   0.006378   3.410 0.000653 ***
## trans.sulphates             0.016897   0.005950   2.840 0.004533 ** 
## trans.alcohol               0.099236   0.007866  12.616  < 2e-16 ***
## trans.colorwhite           -0.003484   0.010721  -0.325 0.745182    
## trans.quality.factorMala   -1.248800   0.013003 -96.036  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.386 on 5175 degrees of freedom
## Multiple R-squared:  0.7401, Adjusted R-squared:  0.7395 
## F-statistic:  1228 on 12 and 5175 DF,  p-value: < 2.2e-16

Analizamos la multicolinealidad mediante el análisis del valor VIF

vif(modelo_2)
##        trans.fixed.acidity     trans.volatile.acidity 
##                   1.277776                   1.228023 
##          trans.citric.acid       trans.residual.sugar 
##                   1.123980                   1.447036 
##            trans.chlorides  trans.free.sulfur.dioxide 
##                   1.444513                   1.864524 
## trans.total.sulfur.dioxide                   trans.pH 
##                   2.356186                   1.242883 
##            trans.sulphates              trans.alcohol 
##                   1.066088                   2.004646 
##                trans.color       trans.quality.factor 
##                   1.000314                   1.299590

Los valores VIF son menores de 5 para todas las variables, en otras palabras la colinealidad es muy baja con el modelo 2.

Ajustamos un modelo nuevo con las variables que presentan significación:

modelo_3 <- lm(trans.quality ~  trans.volatile.acidity +trans.free.sulfur.dioxide +trans.alcohol , new_wine_train) 
summary(modelo_3)
## 
## Call:
## lm(formula = trans.quality ~ trans.volatile.acidity + trans.free.sulfur.dioxide + 
##     trans.alcohol, data = new_wine_train)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -2.43180 -0.43398 -0.01586  0.47378  2.08740 
## 
## Coefficients:
##                            Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                5.816748   0.009104  638.89   <2e-16 ***
## trans.volatile.acidity    -0.140007   0.009994  -14.01   <2e-16 ***
## trans.free.sulfur.dioxide  0.145472   0.010077   14.44   <2e-16 ***
## trans.alcohol              0.382024   0.009762   39.13   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.6538 on 5184 degrees of freedom
## Multiple R-squared:  0.2532, Adjusted R-squared:  0.2527 
## F-statistic: 585.8 on 3 and 5184 DF,  p-value: < 2.2e-16

Los valores VIF son menores de 1.1 para todas las variables, en otras palabras la no hay colinealidad en el modelo 3.

vif(modelo_3)
##    trans.volatile.acidity trans.free.sulfur.dioxide 
##                  1.006974                  1.083543 
##             trans.alcohol 
##                  1.076535

Árbol de decición.

Al igual que para el modelo realizado anteriormente, para la futura evaluación del árbol de decisión, es necesario dividir el conjunto de datos en un conjunto de entrenamiento y un conjunto de prueba.
Esta vez tras la división se llevará a cabo una normalización de los conjuntos de datos para evitar que alguna variable con valores más altos pueda influenciar a las otras.

No hay ninguna proporción fijada con respecto al número relativo de componentes de cada subconjunto, pero la más utilizada acostumbra a ser 2/3 para el conjunto de entrenamiento y 1/3, para el conjunto de prueba.

La variable por la que clasificaremos es el campo de si el vino es de alta o baja calidad, que está en la última columna.

Se divide el conjunto de datos en conjuntos de entrenamiento y prueba. La relación entrenamiento / prueba es normalmente de 70/30 ((2*6486/3=4.324). Queremos asegurarnos de que cada división representa el conjunto de datos en términos de las proporciones de la variable objetivo. La variable objetivo aquí es “Coste”

set.seed(1)
new_wine_arbol<-new_wine[,1:12]
train <- new_wine_arbol[1:4324,]
test <- new_wine_arbol[4325:6486,]
# Filas y Columnas
n = nrow(new_wine_arbol) 
p = ncol(new_wine_arbol)

Normalizamos el set de entrenamiento para que obtener un rango entre 0-1

normalizado_train = function(x) (x - min(x))/(max(x) - min(x))
train.norm = data.frame(apply(train[,-p], 2, normalizado_train), 
                        quality = train[,p])
summary(train.norm)
##  trans.fixed.acidity trans.volatile.acidity trans.citric.acid
##  Min.   :0.0000      Min.   :0.0000         Min.   :0.0000   
##  1st Qu.:0.3649      1st Qu.:0.3640         1st Qu.:0.3617   
##  Median :0.4968      Median :0.4842         Median :0.4681   
##  Mean   :0.5013      Mean   :0.4802         Mean   :0.4888   
##  3rd Qu.:0.6429      3rd Qu.:0.6011         3rd Qu.:0.6170   
##  Max.   :1.0000      Max.   :1.0000         Max.   :1.0000   
##  trans.residual.sugar trans.chlorides  trans.free.sulfur.dioxide
##  Min.   :0.0000       Min.   :0.0000   Min.   :0.0000           
##  1st Qu.:0.2723       1st Qu.:0.3893   1st Qu.:0.3788           
##  Median :0.5773       Median :0.5211   Median :0.5073           
##  Mean   :0.5341       Mean   :0.5014   Mean   :0.5051           
##  3rd Qu.:0.7487       3rd Qu.:0.6241   3rd Qu.:0.6379           
##  Max.   :1.0000       Max.   :1.0000   Max.   :1.0000           
##  trans.total.sulfur.dioxide trans.density       trans.pH     
##  Min.   :0.0000             Min.   :0.0000   Min.   :0.0000  
##  1st Qu.:0.3821             1st Qu.:0.3263   1st Qu.:0.3903  
##  Median :0.4979             Median :0.4571   Median :0.5096  
##  Mean   :0.5073             Mean   :0.4729   Mean   :0.5131  
##  3rd Qu.:0.6365             3rd Qu.:0.6009   3rd Qu.:0.6427  
##  Max.   :1.0000             Max.   :1.0000   Max.   :1.0000  
##  trans.sulphates  trans.alcohol       quality     
##  Min.   :0.0000   Min.   :0.0000   Min.   :4.000  
##  1st Qu.:0.3658   1st Qu.:0.2989   1st Qu.:5.000  
##  Median :0.5052   Median :0.4901   Median :6.000  
##  Mean   :0.4979   Mean   :0.4835   Mean   :5.831  
##  3rd Qu.:0.6204   3rd Qu.:0.6468   3rd Qu.:6.000  
##  Max.   :1.0000   Max.   :1.0000   Max.   :7.000

Normalizamos el set de prueba para que obtener un rango entre 0-1

normalizado_test = function(x) (x - min(x))/(max(x) - min(x))
test.norm = data.frame(apply(test[,-p], 2, normalizado_test), 
                        quality = test[,p])
summary(test.norm)
##  trans.fixed.acidity trans.volatile.acidity trans.citric.acid
##  Min.   :0.0000      Min.   :0.0000         Min.   :0.0000   
##  1st Qu.:0.3509      1st Qu.:0.3640         1st Qu.:0.3404   
##  Median :0.4804      Median :0.4842         Median :0.4255   
##  Mean   :0.4851      Mean   :0.4815         Mean   :0.4557   
##  3rd Qu.:0.6071      3rd Qu.:0.6011         3rd Qu.:0.5532   
##  Max.   :1.0000      Max.   :1.0000         Max.   :1.0000   
##  trans.residual.sugar trans.chlorides  trans.free.sulfur.dioxide
##  Min.   :0.0000       Min.   :0.0000   Min.   :0.0000           
##  1st Qu.:0.2937       1st Qu.:0.3893   1st Qu.:0.3956           
##  Median :0.5980       Median :0.5211   Median :0.5110           
##  Mean   :0.5526       Mean   :0.4951   Mean   :0.5115           
##  3rd Qu.:0.7696       3rd Qu.:0.6241   3rd Qu.:0.6425           
##  Max.   :1.0000       Max.   :1.0000   Max.   :1.0000           
##  trans.total.sulfur.dioxide trans.density       trans.pH     
##  Min.   :0.0000             Min.   :0.0000   Min.   :0.0000  
##  1st Qu.:0.3843             1st Qu.:0.2953   1st Qu.:0.3794  
##  Median :0.4974             Median :0.4334   Median :0.4878  
##  Mean   :0.5110             Mean   :0.4511   Mean   :0.4893  
##  3rd Qu.:0.6401             3rd Qu.:0.5937   3rd Qu.:0.6007  
##  Max.   :1.0000             Max.   :1.0000   Max.   :1.0000  
##  trans.sulphates  trans.alcohol       quality     
##  Min.   :0.0000   Min.   :0.0000   Min.   :4.000  
##  1st Qu.:0.3939   1st Qu.:0.3222   1st Qu.:5.000  
##  Median :0.5135   Median :0.5276   Median :6.000  
##  Mean   :0.5133   Mean   :0.5148   Mean   :5.826  
##  3rd Qu.:0.6306   3rd Qu.:0.6778   3rd Qu.:6.000  
##  Max.   :1.0000   Max.   :1.0000   Max.   :7.000

Construímos el árbol de decisión

arbol = rpart(quality~., data=train.norm)
rpart.plot(arbol)

Modelo modelo no supervisado y basado en el concepto de distancia

Vamos a aplicar un modelo no supervisado. Para consegirlo no usaremos la columna quality.factor, que es la variable que se quiere predecir. Por lo tanto, intentaremos encontrar agrupaciones usando el resto de atributos que caracterizan a cada vino.

new_wine_cluster <- new_wine[,c(1:12)]

Como inicialmente no conocemos el número óptimo de clústers, probamos con varios valores

d <- daisy(new_wine_cluster) 
resultados <- rep(0, 10)
for (i in c(2,3,4,5,6,7,8,9,10))
{
  fit           <- kmeans(new_wine_cluster, i)
  y_cluster     <- fit$cluster
  sk            <- silhouette(y_cluster, d)
  resultados[i] <- mean(sk[,3])
}

Mostramos en un gráfica los valores de las siluetas media de cada prueba para comprobar que número de clústers es el mejor

plot(2:10,resultados[2:10],type="o",col="blue",pch=0,xlab="Numbero de clusters",ylab="Silueta")

Aunque el valor esperado es k=2, dado que el conjunto original tiene 2 factores de calidada (Bueno y Malo), el mejor valor que se obtiene es k=4.

Otro forma de evaluar cual es el mejor número de clústers es considerar el mejor modelo, aquel que ofrece la menor suma de los cuadrados de las distancias de los puntos de cada grupo con respecto a su centro (withinss), con la mayor separación entre centros de grupos (betweenss). Como se puede comprobar es una idea conceptualmente similar a la silueta. Una manera común de hacer la selección del número de clústers consiste en aplicar el método elbow (codo), que no es más que la selección del número de clústers en base a la inspección de la gráfica que se obtiene al iterar con el mismo conjunto de datos para distintos valores del número de clústers. Se seleccionará el valor que se encuentra en el “codo” de la curva

resultados <- rep(0, 10)
for (i in c(2,3,4,5,6,7,8,9,10))
{
  fit           <- kmeans(new_wine_cluster, i)
  resultados[i] <- fit$tot.withinss
}
plot(2:10,resultados[2:10],type="o",col="blue",pch=0,xlab="Número de clusters",ylab="tot.tot.withinss")

En este caso el número óptimo de clústers no está tan claro dado que la curva no llega a estabilizarse.

También se puede usar la función kmeansruns del paquete fpc que ejecutá el algoritmo kmeans con un conjunto de valores y selecciona el valor del número de clústers que mejor funcione de acuerdo a dos criterios la silueta media (“asw”) y Calinski-Harabasz (“ch”).

library(fpc)
#x$duración_corregida<-as.numeric(x$duración_corregida)
fit_ch  <- kmeansruns(new_wine_cluster, krange = 1:10, criterion = "ch") 
fit_asw <- kmeansruns(new_wine_cluster, krange = 1:10, criterion = "asw") 

Podemos comprobar el valor con el que se ha obtenido el mejor resultado y también mostrar el resultado obtenido para todos los valores de k usando ambos criterios

fit_ch$bestk
## [1] 2
fit_asw$bestk
## [1] 2
plot(1:10,fit_ch$crit,type="o",col="blue",pch=0,xlab="Número de clústers",ylab="Criterio Calinski-Harabasz")

plot(1:10,fit_asw$crit,type="o",col="blue",pch=0,xlab="Número de clústers",ylab="Criterio silueta media")

Los resultados son muy parecidos a los que hemos obtenido anteriormente. Con el criterio de la silueta media se obtienen dos clústers y con el Calinski-Harabasz se obtiene un peperfil parecido que muestra 2 clúster.

Como se ha comprobado, conocer el número óptimo de clústers no es un problema fácil. Tampoco lo es la evaluación de los modelos de agregación.

Como en el caso que estudiamos sabemos que los datos pueden ser agrupados en 2 factores de calidad, vamos a ver como se ha comportado el kmeans en el caso de pedirle 2 clústers. Para eso comparamos visualmente los campos dos a dos, con el valor real que sabemos está almacenado en el campo “quality.factor” del dataset original.

new_wine_cluster3 <- kmeans(new_wine_cluster, 2)

par(mfrow = c(1, 2)) 
# trans.residual.sugar y trans.density
plot(new_wine_cluster[c(4,8)], col=new_wine_cluster3$cluster)
# trans.residual.sugar y trans.alcohol
plot(new_wine_cluster[c(4,11)], col=new_wine_cluster3$cluster)

El atributo trans.residual.sugar parece hacer un mejor trabajo para dividir los vinos entre los dos factores de calidad.

2.2.14 Análisis prescriptivo

Regresión Lineal

Predecimos valores utilizando el set de prueba para el modelo lineal.

# Se redondea el resultado para eliminar los decimales en la predicción
modelo_3_Pred = round(predict(modelo_3, newdata=new_wine_test))
head(modelo_3_Pred)
## 12 16 25 34 58 61 
##  5  6  5  6  5  5

Creamos la matriz de confusión para comprobar la precisión del modelo. Se enfrentan los valores predichos por el modelo con los valores reales del set de prueba.

confusionMatrix(table(modelo_3_Pred, new_wine_test$trans.quality))
## Confusion Matrix and Statistics
## 
##              
## modelo_3_Pred   4   5   6   7
##             4   2   0   0   0
##             5  21 159 100  10
##             6  19 216 502 243
##             7   0   0  11  15
## 
## Overall Statistics
##                                           
##                Accuracy : 0.5223          
##                  95% CI : (0.4948, 0.5498)
##     No Information Rate : 0.4723          
##     P-Value [Acc > NIR] : 0.00017         
##                                           
##                   Kappa : 0.1689          
##                                           
##  Mcnemar's Test P-Value : NA              
## 
## Statistics by Class:
## 
##                      Class: 4 Class: 5 Class: 6 Class: 7
## Sensitivity          0.047619   0.4240   0.8189  0.05597
## Specificity          1.000000   0.8581   0.3022  0.98932
## Pos Pred Value       1.000000   0.5483   0.5122  0.57692
## Neg Pred Value       0.969136   0.7857   0.6509  0.80110
## Prevalence           0.032357   0.2889   0.4723  0.20647
## Detection Rate       0.001541   0.1225   0.3867  0.01156
## Detection Prevalence 0.001541   0.2234   0.7550  0.02003
## Balanced Accuracy    0.523810   0.6410   0.5606  0.52265

Nuestro modelo tiene una precisión del 52.23%

Diagnosticamos el modelo a continuación

#Revisamos los residuos para asegurar que están distribuidos normalmente
hist(modelo_3$residuals)

#Comprobamos la media de los residuos para asegurarnos que ésta se encuentra cerca de 0
mean(modelo_3$residuals)
## [1] 1.041235e-17
#Visualizamos los diagramas
layout(matrix(c(1,2,3,4),2,2))
plot(modelo_3)

Tras examinar los resultados podemos afirmar: * Se satisface la condición de linealidad.
* El análisis gráfico confirman la normalidad.
* No hay evidencias de falta de homocedasticidad.
* No hay multicolinialidad.

Árbol de decisión

Predecimos valores utilizando el set de prueba para el modelo del árbol de decisión.

predicted_model <- round(predict( arbol, test.norm[,-p] ))
head(predicted_model)
## 6527 6528 6529 6532 6533 6535 
##    6    6    5    6    5    6

Una vez tenemos el modelo, podemos comprobar su calidad prediciendo la calidad para los datos de prueba que nos hemos reservado al principio.

print(sprintf("La precisión del árbol es: %.4f %%",100*sum(predicted_model == test.norm) / length(predicted_model)))
## [1] "La precisión del árbol es: 56.1517 %"

La precisión del modelo del árbol de decisión es semejante al modelo lineal anterior.

Clusterización

Tras crear el modelo se procede a comprobar cuales de los vinos han sido mal clasificados y cómo. Eso lo conseguimos con el siguiente comando.

table(new_wine_cluster3$cluster,new_wine$trans.quality.factor)
##    
##     Buena Mala
##   1  2898  802
##   2  1464 1322

Los resultafos de la tabla muestran que el cluster 1 corresppnde a Bueno y el dos a Malo. Y así, podemos sacar un porcentaje de precisión del modelo:

Bien calsificados 2900 + 1320 = 4.220 Mal clasificados 1462 + 804 = 2.266

print(sprintf("La precisión de la clusterización es: %.4f %%",100*4220/(4220+2266)))
## [1] "La precisión de la clusterización es: 65.0632 %"

3 Representación de los resultados a partir de tablas y gráficas

Durante el análisis exploratorio de los datos se han mostrado múltiples representaciones, de entre todas ellas mostramos a modo de reumen la relación entre las variables.

library(psych)
pairs.panels(new_wine[,c(1:12)])

Los resultados muestran simetría en todas los atributos excepto en dos. La simetría se consiguió aplicando la transformación de Boxcox para normalizar todas las variables. Las dos únicas variables que no muestran simetría son trans.citric.acid y trans.density y son las dos únicas que hemos representado mediante diagramas Q-Q para observar su distribución. Representar el resto de variables transformads no hubiera mostrado ningún resultado particularmente interesante, ya que están normalizadas y lo comprobamos en la diagonal de la representación que mostramos en este apartado.

4 Conclusión

Tras realizar la limpieza de los datos y aplicar distintos modelos para averiguar la calidad de los vinos concluímos que el modelo de clusterización es el que más precisión muestra. Sin embargo, hay que resaltar que para la construcción de los modelos no se ha llevado una procedimiento uniforme, ya que con el objetivo de practicar distintos enfoques, se han usado distintas técnicas para normalizar, estandarizar, homogenizar y separar los datos en test y prueba.

Es evidente, que para poder tomar una conclusión acertada respecto al los resultados obtendos debería haberse seguido un proceso más uniforme, pero como se ha comentado, el objetivo es practicar.

5 Archivo preprocesado

Para finalizar se genera un archivo con el conjunto de datos procesado para análisis posteriores.

# Guardamos los datos
write.csv(new_wine, file="Wine_clean.csv",row.names=F)

6 Repositorio

A continuación se muestra la dirección de acceso al repositorio, donde podemos encontrar el código fuente, el archivo original, el archivo preprocesado y el documento html (RMarkdown).

https://github.com/SevillaFe/DataCleaning_DataAnalysis_UOC

7 Tabla de contribuciones y referencias

Tanto la la tabla de contribuciones como las referncias utilizadass durante esta práctica, se encuentran disponibles en la wiki del repositorio.